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While a rich variety of self-propelled particle models propose to explain the collective motion of 
fish and other animals, rigorous statistical comparison between models and data remains a chal¬ 
lenge. Plausible models should be flexible enough to capture changes in the collective behaviour of 
animal groups at their different developmental stages and group sizes. Here we analyse the statis¬ 
tical properties of schooling fish {Pseudomugil signifer) through a combination of experiments and 
simulations. We make novel use of a Boltzmann inversion method, usually applied in molecular dy¬ 
namics, to identify the effective potential of the mean force of fish interactions. Specifically, we show 
that larger fish have a larger repulsion zone, but stronger attraction, resulting in greater alignment 
in their collective motion. We model the collective dynamics of schools using a self-propelled par¬ 
ticle model, modified to include varying particle speed and a local repulsion rule. We demonstrate 
that the statistical properties of the fish schools are reproduced by our model, thereby capturing a 
number of features of the behaviour and development of schooling fish. 

PACS numbers: 87.18.-h, 05.65.+b, 05.10.-a, 05.40.-a 


INTRODUCTION 


In sufficiently large collective systems, the behaviour 
of an individual can be dominated by the generic statis¬ 
tical effects of many individuals interacting, rather than 
its own behaviour [l[. Much of the progress in under¬ 
standing collective motion of animal groups has involved 
applying ideas borrowed from the statistical physics of 
materials like magnets or fiuids For example, 

changes in group densities produce phase transitions at 
critical group sizes [3,[1]- More complex collective states, 
such as swarm, mills and polarised groups depend on the 
density of a group and the noise within the system [^. 
Recent studies of starlings and midges have looked at 
spatial velocity fiuctuations EEl, long range correla¬ 
tions fl^. an d diffusive [l^ and entropic characteristics 
of hocks [l^. Other experiments with artificial particles 
have looked for similarities and differences between self- 
organised living matter and thermal equilibrium systems 
[l5j. These latter approaches gather statistical informa¬ 
tion about self-organising structures in order to param- 
eterise models (see for example the maximum entropy 
approach [3, [I^). However, none of these have explic¬ 
itly solved the inverse problem of using the macro-level 
properties of animal groups to find out how the individ¬ 
uals within them interact. 


ter and many of its properties 17|. At the same time, 
molecular interaction potentials are difficult to measure 
experimentally and hard to compute from first prin¬ 
ciples. An alternative approach, therefore, is to esti¬ 
mate them from experimentally determined structures of 
molecules. The interactions in these structures are usu¬ 
ally strongly coupled and assemblies are typically driven 
by weak forces (e.g., hydrophobicity or entropy) 0. 
Therefore, estimation of these potentials requires appli¬ 
cation of sophisticated coarse-grained techniques such as 
reverse Monte-Carlo [l^ , inverse Monte Carlo or It¬ 
erative Boltzmann inversion [El- These methods adjust 
the force field iteratively, until the distribution functions 
of the reference system are reproduced as accurately as 
possible. In other cases, when the potentials are uncou¬ 
pled or weakly coupled, a more straightforward direct 
Boltzmann inversion approach can be applied, which ap¬ 
proximates the potential by the negative logarithm of the 
radial distribution function (Ej |. In collectively moving 
animal groups the interactions between members are usu¬ 
ally assumed to be of hierarchical structure, with repul¬ 
sion having highest priority at small distances . Thus 
one can expect that the latter method can be also applied 
to animal self-organised systems, such as fish schools, to 
infer the interactions within these groups from experi¬ 
mental data. 


This inference problem is essentially a statistical 
physics problem. The last few decades have seen a major 
increase in research at the interface of molecular dynam¬ 
ics and biophysics. In soft matter systems, estimating the 
potential energy of an interaction and the corresponding 
potentials is of particular importance, as the strength of 
intermolecular interactions determines the state of mat¬ 


Here we investigate the schooling behaviour of fish us¬ 
ing Boltzmann inversion and related methods. Unlike 
molecules and physical particles, fish change their be¬ 
haviour as they go through various developmental stages 
[ 2 ^ . For example, onset of schooling is only possible 
when the central nervous system of fish is sufficiently de¬ 
veloped to support a high level of coordination of visual 
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and mechanosensory information 23- The developmen¬ 
tal differences are usually also reflected in changes of the 
key characteristics of motion, including speed. There¬ 
fore, fish of different sizes can not be considered simply as 
particles of different physical size, since their behaviour 
changes with their size. We thus expect the statistical 
properties of the group, and of individuals, to change 
both with the density of fish and their developmental 
stage. 


METHODS 
Experimental details 

We used groups of 10 to 60 Pacific blue-eyes {Pseudo- 
mugil signifer) with approximately three different body 
lengths (from hereafter referred to as fish sizes): ~ 7.5 
mm (small), ^ 13 mm (medium) and ~ 23 mm (large) 
(see the electronic supplementary material figure SI) . 
Because body size is related to the age of a fish [26j, 
the three body lengths used in this study likely repre¬ 
sent three distinct age classes. The largest fish (23 mm) 
constituted of sexually mature individuals, although we 
observed no sexual behaviour in the trials. The fish were 
confined into a large shallow circular arena (760 mm di¬ 
ameter) and filmed from above at high spatial and tem¬ 
poral resolution. The positions of fish were subsequently 
tracked using DIDSON tracking software [l|. On aver¬ 
age 86 % of fish were identified and tracked in our exper¬ 
iments, which is a similar level of accuracy as compared 
to other studies that track large numbers of individuals 

i- 


Model 


In our two-dimensional model, the fish are represented 
by N point particles at number density p and variable 
particle speed Vi. The system undergoes discrete-time 
dynamics with a time step At. The direction of motion 
of each particle (Fig. [J) is affected by repulsive or align¬ 
ment interactions with other particles located inside the 
zone of repulsion (zor) or zone of alignment (zoa), respec¬ 
tively. Time evolution therefore consists of two steps: ve¬ 
locity updating and streaming (position update). In the 
first computational step, position of each particle (r^) 
is compared to the location of the nearest neighbours. 
The repulsion rule has an absolute priority in the model 
and is modelled as a typical collision avoidance [H, 1 ^ 


u{t)i 


Tir 


-E i'd(t)/ 


fLj' 




with Vij = Tj — Ti and Ur 


being a number of particles inside zor. The alig nment 
rule similar to one used in the Vicsek model (29j| takes 
into account velocities of all particles located inside the 



FIG. 1. The illustration of the interaction parameters in the 
SPP model. The particle shown in red turns away from the 
nearest neighbours within the zone of repulsion (zor) to avoid 
collisions and aligns itself with the neighbours within the zone 
of alignment (zoa). 


zone of alignment u{t)i = ^ Vj(t)/ ^ Vj(t) with 

1=1 1=1 

being a number of particles inside zoa. The velocities of 
particles are updated according to 


V,(t) = z;,(t)u(t),i?i(e.(t))i?2(0*(t)) (S7) 

with Vi{t) = vo[tp(t)]'^ defining the particle individual 
speed Vi{t) based on the averaged local order in¬ 
side both behavioural zones [30|-l33| . Vi (t) takes its maxi¬ 
mal value Vi{t) = vq when velocities of particles inside 
zor and zoa are perfectly aligned = 1 while ab¬ 

sence of local order '0(t) = 0 results in Vi{t) = 0. The 
exponent 7 controls the sharpness of the speed change. 
Note that for any 7 an isolated particle will move with 
maximal speed vq. The misaligning noise is introduced 
through a random rotation of the resulting 

particle velocity according to a Gaussian distribution 

where ^i(t) is a random 
variable and rj is the noise strength. 

Wall avoidance is modelled as a particle orientation 
adjustment through rotation R 2 { 0 i{t)) of the particle 
velocity with a time-dependent turning rate 0i(t) = 
vo(j>i{t)/di{t) . is the angle between the heading of a 

fish and normal to a time-dependent point of impact on 
the wall [33, IHj . di (t) denotes a distance from particle 
i to the impact point. Such construction of the rotating 
rate allows to achieve its strong damping at large dis¬ 
tances from the wall and for smaller angles of approach 
of a collision point on the wall. At these conditions its 
influence on particle’s motion is insignificant. 

When the velocity update step is complete, the particle 
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FIG. 2. Example time evolution of the spatial distribution in 
arena at = 10 small fish (a.), N — 60 small fish (b), and 
= 10 large fish (c). Each cell in plots (a)-(c) denotes a 
20 degrees radial segment of the arena and 1 second of time. 
0 is the angle measured counterclockwise from the positive 
direction of the a;-axis as defined by the camera position. 


positions are updated by 

ri{t +At) =n{t)+Y^{t)At (1) 

The model was parametrised based on the experimen¬ 
tally obtained data. The unit of length in our simula¬ 
tions is equivalent to the metric length used in the ex¬ 
periment. To set the unit of time we choose a particle 
speed vq = bv^, where b is the behavioural reaction time 
[3q of fish (6 = 0.05 s) and Vq is the average speed in 
experiment. The integration time step was set to At = 1 
for all simulations. The noise strength rj was fixed at 0.1 
for all trials. The exponent 7 was set to 1 for in all sim¬ 
ulations. Fish of different size are modelled by scaling 
the size of the alignment zone with a factor k propor¬ 
tionally to experimentally measured differences in body 
lengths so that k = 1, k = 1.73 and k = 3.07 correspond 
to small (7.5 mm), medium (13 mm) and large (23 mm) 
fish, respectively. 

Total number of time steps in each run was 1 x 10®. 
The statistics was collected in the steady state and each 
characteristic of motion was calculated by averaging over 
five independent runs. The radius of the arena was fixed 
at i? = 380 for all simulations. The initial conditions for 
fish positions and velocities were chosen at random from 
the uniform distribution. 


RESULTS 

We first investigated the spatial distribution of fish in 
the arena. Small fish displayed limited collective motion. 
For example, 10 small fish tended to form a dispersed 
group, where most of the fish moved very little (figure 
[2Ka)). Larger groups of 60 small fish showed slightly more 
collective motion, but not all fish moved in the same di¬ 
rection at the same time (figure [2(b)). In contrast, even 


small groups of large fish showed highly aligned collective 
motion (figure | 2 (c)). 

We calculated the average area. A, covered by the 
group using a convex hull algorithm (see the electronic 
supplementary material). The average value of A is plot¬ 
ted for different group and fish sizes (figure |3(a)). For all 
three fish sizes, larger groups occupied a larger area. The 
density, p = N/A, also increased with the number of fish 
in the group (figure | 2 (b)), suggesting that the fish pack 
closer together in larger groups. Figure |3(a) indicates 
that groups of small fish occupied a larger area than the 
groups of medium-size or big fish. This finding supports 
the results of the spatial analysis (figure [ 2 ]); small fish 
were more dispersed over the arena. As a result, groups 
consisting of small fish were less dense (figure |3(b)) than 
groups of larger fish. 

To better quantify the spatial arrangement of groups, 
we measured their packing fraction a (figure |3(c)). This 
is the ratio between the total body area of all fish in a 
group {Af = the global area of a group 

(A): a = AflA, where Ai is the body area of individual 
fish. For all body sizes of fish, packing fraction increased 
with group size. Groups of smaller fish had the low¬ 
est packing fraction ranging between 0.001 and 0.004 for 
groups of 10 and 60 individuals respectively. In contrast, 
groups of medium-size and large fish had higher packing 
fractions of a > 0.043 and a > 0.054 respectively. The 
lowest packing fractions in groups of small fish are com¬ 
parable to those observed in bird flocks [H, , whilst 

the larger packing fractions approach those of some bac¬ 
teria (39|. In physical systems, small values of a typi¬ 
cally correspond to gases, while larger values (a > 0.4) 
to liquids or crystals [i^. All packing fractions observed 
in our experiments, therefore, are comparable with an 
atomistic system in its gaseous state. At the same time, 
the large differences between packing fractions for small 
and medium-size fish and for small and large fish, reach¬ 
ing one order of magnitude (t-test, p < 1 x 10 “®), suggest 
possible differences in other statistical characteristics of 
the system for varying fish size. 

We next characterised ordering in our system using the 
polar order parameter 


if = 



N 


exp(z6»i) 


( 2 ) 


where i is the imaginary unit, 9i is the direction of mo¬ 
tion of individual fish, and (•) denotes the time average. 
We should note that the polar order parameter is gener¬ 
ally sensitive to the choice of confining geometry within 
which the agents move. For instance, for a small circular 
arena, the relationship between the radius of the arena 
and radial positions and speeds of the agents determines 
the maximum value of polarisation that can be reached 
in the system. In our experiments, we used a large arena 
with a radius exceeding the average group width for large 
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FIG. 3. Statistical properties for groups of fish with three 
different body lengths: 7.5, 13, and 23 mm. (a) Area of a 
group A. (b) Number density p. (c) Packing fraction a. (d) 
Polar order parameter ip. 

fish by more than 3 times. Large fish also occupied on av¬ 
erage a radial segment of only 30 degrees (see figure[2Kc)). 
Thus we expect that the polar order parameter can pro¬ 
vide meaningful results in characterisation of alignment 
in our experimental system. Values of ip are plotted in 
figure [SKd) as a function of group and fish size. Small 
fish (7.5mm) did not display much ordering for any group 
size with p k, 0.4 — 0.55, confirming previous observations 
(figures [2Ia),(b)). Groups of medium-size fish were the 
most ordered, with values of the order parameter ranging 
between 0.73 and 0.84 depending on group size. Large 
fish displayed slightly lower order than medium-size fish. 

We then investigated the nature of the interactions be¬ 
tween the fish using statistical mechanics. We started 
by looking at the pair distribution function (PDF) 
which allowed us to study how the local density around 
each fish varied with respect to the average density in the 
system. It is defined by 

IN{r)N{r) \ 

where 6 is the Dirac delta function, |ry| is the distance 
between fish i and j, S{r) is the surface area of a shell, 
N{r) is the number of fish inside a shell and (•) stands for 
the time average (see the electronic supplementary ma¬ 
terial for details). A set of PDF-curves g{r) for various 
fish sizes is presented in figure SKa). Small fish tended 
to form aggregations with densities up to 4 times above 
the average density of the system and with a maximal 
half-radial width of more than 25 fish body lengths. For 
medium and large fish, the maximum density in a clus¬ 
ter was 8 times larger than the average in the system. 
The size of the aggregation of medium and big fish was 
as large as 17 and 10 body lengths, respectively. An¬ 
other notable difference between the three curves is the 
location of the local density peak. For small fish, the 


Experiments 



r, mm r, mm r, mm 


FIG. 4. Statistical properties of fish in experiments, (a) Pair 
distribution function g{r) smoothed with 5-point moving av¬ 
erage, (b) effective potential of the mean force of the inter¬ 
action U{r), and (c) mean force of the interaction F{r) for 
small (7.5 mm) medium (13 mm), and large (23 mm) fish, 
(d) Pair distribution function g{r), (e) effective potential of 
the mean force of the interaction U{r), and (f) mean force of 
the interaction F{r) for medium-sized fish (13 mm) in groups 
of different size. The top row contains plots for different body 
size and the bottom row contains plots for different group size. 
The top and bottom legends correspond to plots in the top 
and bottom rows, respectively. 

peak is at 24.5 mm, whereas for larger fish it is signif¬ 
icantly shifted towards 30 - 40 mm. Figure lU^d) shows 
PDF plots for medium-sized fish at varying group size. 
For the smallest groups of 10 fish, the maximum density 
observed is 25 times above the system’s average density. 
The peak value of g{r) decreases with increasing group 
size (from 25 for groups of 10 fish down to 6.5 for the 
largest groups of 60 fish) whilst the position of the maxi¬ 
mum remains unchanged at approximately 30 mm for all 
group sizes. The maximal half-radial width of the aggre¬ 
gation increases with group size from 145 mm for groups 
of 10 fish to 240 mm for groups of 60 fish. 

The differences in the pair distribution function sug¬ 
gest there is large variation in the underlying pair po¬ 
tential. In our system, in the absence of any external 
field, this potential is the effective potential of the mean 
force of the interactions between fish. Studies of active 
matter show that the resulting steady states of such sys¬ 
tems often satisfy Boltzmann distribution , even given 
that these active systems are essentially out of equilib¬ 
rium. Here we apply the opposite route: we start with 
the assumption that the steady-state configurations ob¬ 
served in the experiment are drawn from the Boltzmann 
distribution 

P{r) = Z~^ exp[—pU{r)], (4) 

where P = l/ksT = 1 with a Boltzmann constant 
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FIG. 5. Statistical properties of a model system, (a) Pair 
distribution function g(r), (b) effective potential of the mean 
force of the interaction U(r), and (c) mean force of the in¬ 
teraction F{r) for small (k = 1), medium (k = 1.73), and 
large (fc = 3.07) simulated fish, (d) Pair distribution function 
g(r), (e) effective potential of the mean force of the inter¬ 
action U{r), and (f) mean force of the interaction F{r) for 
medium-sized fish (13 mm) in groups of different size. The 
top row contains plots for different body size and the bottom 
row contains plots for different group size. The top and bot¬ 
tom legends correspond to plots in the top and bottom rows, 
respectively. 


kB and the system’s temperature denoted by T; Z = 
f exp[—/3f7(r)]dr is the partition function and U(r) is 
the effective potential of the interaction. In this work, 
the choice of the inverse temperature /3 = 1 is arbitrary 
as this term is not related to thermodynamic tempera¬ 
ture in our system. Instead, it accounts for temperature¬ 
like fluctuations in the system which we assume to take 
the same value for all three size classes of fish considered 
here. The fact that we observe large differences in the 
stationary distributions of different sized fish and differ¬ 
ent group sizes is related to the changes in motion of in¬ 
dividuals when the configuration of neighbours changes. 
Nevertheless, the average number of individuals within 
each shell of the pair distribution function remains con¬ 
stant over the whole duration of an experiment, defining 
a steady-state in our system. 

To derive the effective potential of the mean force of 
the interactions we use the direct Boltzmann inversion 



[/(r) = —kBTlng(r). (5) 

Figures |4Kb),(e) present the effective potential energy 
profiles for different sized fish and different group sizes, 
respectively. Note that all the curves on both figures have 
practically the same slope for the decreasing part of the 
potential down to I7(r) = —0.55. The minimum of these 
curves occurs at a greater distance for larger fish (figure 


HKb)), indicating that larger fish have a larger repulsion 
zone. The increasing portions of the curves have sim¬ 
ilar slopes only for medium-sized and large fish (figure 
SKb)). This suggests that fish of different body sizes and 
in groups of various size have similar repulsion strength 
(or collision avoidance potential) at short distances, but 
different attraction strengths towards neighbours. To get 
a conclusive picture of the variation of the interaction 
strength over the separation distance between the indi¬ 
viduals, we calculate the mean force of the interaction 
Fir): 




( 6 ) 


In differently-sized fish (figure IHc)), the repulsive force 
(positive and short portion of the Fir) curves) is stronger 
than the attractive one (negative and long portion of the 
Fir) curves). In other words, repulsion is independent of 
body size and spans a much shorter distance than attrac¬ 
tion. Constant attraction force at large distances in our 
system arises because fish can perceive conspecifics over 
the entire arena. In other systems, a change of the at¬ 
traction force over distance could be used for identifying 
topological interactions between individuals. Previous 
studies ^ used specific channel confinement that mim¬ 
ics a geometrically frustrated anti-ferromagnet to distin¬ 
guish between different types of interactions. Our ap¬ 
proach is simpler and is not limited to specific experi¬ 
mental setups. For example, an abrupt change in the 
strength of attractive force would correspond to metric- 
type interactions. 

To validate the fish interactions established from the 
experiments, we simulated the collective motion of the 
fish using a two-dimensional metric self-propelled parti¬ 
cle model that accounts for variable fish speed and geo¬ 
metrical confinement (see Methods for model description 
and the electronic supplementary material for details on 
other motion statistics). Figures [5l(a)-(c) present plots 
of PDF, effective interaction potential and mean force of 
the interaction for the three different sizes of simulated 
fish. Overall these qualitatively match the experimental 
data. Medium (fc = 1.73) and large (fc = 3.07) simulated 
fish form more dense groups than the small simulated fish 
(fc = 1) with a density of 12.5, 14 and 10 times above the 
average, respectively (figure [SKa)). The density peak is 
also observed at larger distances for bigger simulated fish. 
The repulsive force is practically the same for all three 
cases at any given distance (figure EKc)) as the radii of 
the repulsion zone are constant for all fish sizes. The at¬ 
tractive portion of the F'(r)-curves has a complex shape 
as in the experiment indicating strongest attraction at 
r Ki 28, r « 90 and r Ki 200 mm for small, medium 
and large simulated fish, respectively. For all five group 
sizes of the medium sized simulated fish (fc = 1.73) the 
maximum of the local density is well above the system’s 
average (figure EKd)) and is largest for groups of 10 fish 
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{g{r) = 52). All the curves cross the g{r) = 1 line in the 
same order as in the experiment: The average half-radial 
width of groups of 10 and 60 simulated fish are ^100 mm 
and 210 mm, respectively. As in the experiment, the 
repulsive force (figure EDf)) decays with a distance and 
takes very similar values for all five cases. The attraction 
in groups of 10 and 20 individuals is maximal at r = 90 
and r = 115, respectively, then decreases steeply with 
a distance and reaches zero at r r; 190. For the other 
three curves (N=30, 40 and 60), the attractive force has 
a maximum at r Ri 100. At all intermediate distances the 
dissipation is stronger in small groups of simulated fish 
which is in agreement with experiments. 


DISCUSSION 


Although many recent studies have focused on under¬ 
standing interactions between fish in schools ii 
a detailed description of large systems consisting of tens 
or hundreds of individuals thus far remained a challenge. 
In analysis of such systems, the quality of fish tracking 
becomes a limiting factor, as many of the techniques tra¬ 
ditionally used to study fish behaviour require high con¬ 
sistency of individual fish identities over time in order 
to reconstruct velocity or acceleration profiles of the in¬ 
dividuals. When the long-time individual identities are 
not available, an approach that relies only on individuals’ 
positional data, such as Boltzmann inversion, can be use¬ 
ful for assessing interactions in large groups. The latter 
approach is also faster as it works with spatial data only, 
and the pair distribution function used in its first step 
allows also to extract useful information about aggrega¬ 
tions and clustering in the system. Such approach can 
also potentially simplify fish tracking and post-processing 
of experimental data. 

This approach, however, also has its limitations. The 
direct Boltzmann inversion has become a popular method 
to derive interaction potentials across different fields first 
of all due to its simplicity and general applicability. Even 
though this method has a straightforward nature, the po¬ 
tentials derived with the Boltzmann inversion are gen¬ 
erally non-unique and can be state-dependent. They 
can also be influenced by long-range correlations or by 
the anisotropic nature of the interactions that is inher¬ 
ent to many living systems. Moreover, these potentials 
are effective as they may include multi-particle correla¬ 
tions (e.g. simultaneous response to positioning or move¬ 
ments of several conspecihcs) and correlated contribu¬ 
tions from the surroundings, such as geometrical confine¬ 
ment effects. Such potentials thus do not share the typi¬ 
cal properties attached to equilibrium potentials and this 
limits their use mostly to qualitative description of the 
interactions. This also suggests that to draw a conclu¬ 
sive quantitative picture on how individuals interact in a 
group, a method like Boltzmann inversion needs to be ac¬ 


companied by rigorous analyses of other motion statistics 
and/or by extensive computer simulations of the system 
of interest. 

The challenges of using self-propelled particle models 
for modeling collective animal behaviour are also known 
[ 4 ^ IS ]. The tractable models are usually oversimpli¬ 
fied and too general and thus lack flexibility required to 
reflect the features of a particular phenomena. In this 
work, we tried to overcome many of these limitations by 
parameterising our simple SPP model with experimen¬ 
tal data. Some of the simplifications, on the other hand, 
were introduced in the model deliberately, such as the 
absence of an explicit attraction rule for the interaction 
between the agents (see also Model section). We found 
the inclusion of the attraction rule unnecessary since a 
combination of fish-fish alignment and fish-wall interac¬ 
tions proved to be sufhcient to reproduce the dynamics 
observed in experiments both statistically and visually 
[ 4 ^ . We should note, however, that the effective cohe¬ 
sion detected in the model system by the direct Boltz¬ 
mann inversion method is hence attributed to a combi¬ 
nation of the inter-individual and agent-wall interactions. 
Even more detailed description of the system and possi¬ 
bly a better fit to experimental data could be achieved 
if mass and inertial forces of the individuals are taken 
into account. These features can be relatively easily im¬ 
plemented in a model which describes the motion of in¬ 
dividuals by including friction and stochastic forces, e.g. 
Active Brownian Particle model with aligning interac¬ 
tions [ 45 I. l46l|. 


CONCLUSION 


Using a combination of Boltzmann inversion and tra¬ 
ditional statistical methods, we have inferred how fish 
{Pseudomugil signifer) interact in schools. Previous stud¬ 
ies have applied a force-matching approach to infer the 
interactions of schooling fish from their movements ( 4 ^ . 
Our method can infer these interactions directly from 
the static spatial distribution of individuals in groups. 
Whilst repulsion forces had the same strength for dif¬ 
ferent sized fish, attraction strength increased in larger 
fish, consistent with how a fish’s movement develops with 
age 23 ■ The interactions between fish also changed as a 
function of group size, as suggested by other studies [s^ . 
Our model, refined on the basis of these observations, 
could capture the dynamics of schooling fish. These find¬ 
ings are also in line with the results of our previous study 
44| . where we utilised an observational test to cross- 
validate the model used in the present study. We expect 
that our findings could also generalise to other species 
that exhibit schooling behaviour [13, . Application of 

the approaches used in statistical physics, coupled with 
informed models of collective motion, now allows us to 
shed more light on the intricacies of how individuals in 
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groups interact. 
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EXPERIMENTAL DETAILS 


Materials 


Pacific blue-eye fish {Pseudomugil signifer) were caught in hand nets from Narrabeen Lagoon, New South Wales, 
Australia (33°43'03 S, 151°16'17 E). Fish were kept in filtered freshwater in 150 1 glass tanks at 22 — 25° and fed 
crushed flake food ad libitum. All fish were housed for at least three weeks prior to experimentation. The experimental 
arena was circular with a diameter of 760 mm. It was filled to a depth of 70 mm with aged and conditioned tap 
water. The arena was lit by fluorescent lamps and was visually isolated. For each trial, we randomly selected N fish 
{N = 10,20,30,40,50 or 60 for small and medium fish and N = 10,20,30,40 for large fish) of similar size (see Fig. 
EH) from the housing tanks and placed them in the experimental arena. Fish were left to acclimate to the arena for 
at least five minutes, after which they were filmed for 15-20 minutes at 15 frames per second. We used a Logitech 
Pro 9000 camera placed orthogonally to the arena above its geometrical centre at a distance > 1 m minimising the 
radial distortion. The number of trials for each group size ranged between 3-10 (see Table |T| due to limitations in 
the number of fish we could obtain for large or small body sizes. Because of the large numbers of fish we used for 
the experiment, we reused fish between trials. Fish were never used more than once per day and fish were used a 
maximum of 5 times. 
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FIG. SI. Average fish body length ±1 s.d. for each size class (left). Average fish body area ±1 s.d. for each size class (right). 


Data collection and acquisition 

Films were recorded in .mov format using original camera manufacturer software and subsequently converted to 
.avi using DirectShowSource and VirtualDub (v 1.9.2). The tracking was performed using DIDSON tracking program 
M- The raw data consisted of x and y coordinates, fish identity and a time stamp. The accuracy of the tracking 
process was checked by projecting the raw tracking data onto experimental videos. 

Figure E2] shows radial distribution of fish in the arena calculated as g{R) = p{R)/{p), where p{R) is the density of 
fish in the circular shell of mean radius R and (p) is the average fish density in the system. Small fish I = 7.5 mm are 
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Group size 

Average body length 

7.5 mm 

13 mm 

23 mm 
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40 
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50 
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60 
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TABLE I. Number of videos recorded for each group size and body size. Each film is 15-20 min long. 



FIG. S2. Distribution of fish across the arena in radial direction {N = 40). 0 mm corresponds to the center of the arena; 
the wall is located at 380 mm from the center. The dotted black line at g{R) = 1 corresponds to completely homogeneous 
distribution of fish throughout the arena. 


distributed more regularly throughout the tank as compared to large individuals I = 23 mm tending to move close to 
the arena wall. 


MOTION STATISTICS ANALYSIS 


Calculation of the shell area for the pair distribution function in confinement 


When calculating the pair distribution function for confined systems particular care needs to be taken when the 
particles are located close to the wall, at a distance smaller than the radius of the largest shell. In experiments with 
fish this scenario is common (see Fig. [S2]for illustration). In such situations only area of the shell lying inside the 
boundaries of the confining geometry should be considered. For all cases the shell area can be calculated as a difference 
between the areas of two neighbouring circles. The two latter areas for convenience can be computed separately as 
intersections excluding the areas outside the constraint. If confinement is circular (see Fig. IS3|) . as is in our case, to 
find the intersection area we can use a formula for the circular segment of triangular height d' (excluding the height 
of the arced portion) and radius R' : 


SiR',d') = R'^cos-\d'/R') - d'^jR'^ - (SI) 

In the simplified case that we consider here (with only one circle for the pair distribution function) we have two 
such segments, one for the circular confinement (experimental arena) with a radius R and hight di, and another one 
for the shell boundary of the pair distribution function, having radius r and height ^ 2 . The heights are calculated as 

d^ — + R^ 

di = 


2d 


(S2) 
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FIG. S3. Illustration of the intersection of the pair distribution function circle r with confining circle R. Only one circle for 
the pair distribution function is drawn for simplicity. 


and 


do = 




2d 


(S3) 


To calculate the total area of the intersection Eq. (ISII) needs to be solved two times, once for each segment. Thus, 
combining Eqs. (IS1IIS3I) we get 


S{r') = S{R,di) + S{r, dz) = 
" d^ + 


= cos ^ 


2dr 


+ B? cos 


-1 


d^ + 

2Zr 


— -\/{—d + r + R){d + r — R){d — r + R){d + r + R). 

(S4) 


Eigure[S4] shows a plot of the pair distribution of par¬ 
ticles for highly homogeneous system in circular confine¬ 
ment (distribution of particles is uniform). The black 
curve displays a clear linear decay of g{r) with increasing 
inter-particle separation distance r. This represents the 
case where the area outside the constraint has also been 
included in calculations. The red curve represents the 
case when all the shell areas have been calculated with 
the method described above. This shows highly regular 
distribution g{r) of particles for all separation distances r 
(1:1 relation between the local shell density and average 
density in the system). 


Calculation of the surface area of a fish group and 
the body area of an individual fish 


The surface area of a group was calculated for every 
frame of a video. Eirst, we computed a convex hull C of 
a set of N points representing geometrical centres of fish 
bodies S3| . The convex hull is defined by 


N N 

C = ^ XiPi : A > 0 for all i and ^ Xi = 1, (S5) 

i=l i=l 


where i is the point (fish) index of a point pi and A is 
the non-negative weight coefficient. The resulting convex 
hull gives the identities of the vertices of a polygon as an 
output. The area A of this polygon can be computed as 
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FIG. S4. Pair distribution function g{r) for a homogeneous 
test system constituting of 1000 particles (average over 5 x 10® 
positional configurations), R = 380. 



FIG. S5. Complex shapes obtained from the experimental 
data. Left: Typical polygon returned by a convex hull algo¬ 
rithm. Right: Fish body shape estimated as a convex shape. 
The displayed frame is taken from a video for 40 medium- 
size fish. The frame on the left has been scaled for clarity 
purposes. 


S4| 

1 ^ 

-4 = - ^ (xiy^+i - x^+iyi), (S6) 

^ i=l 

where xt and yi represent coordinates of the vertices of 
a polygon. The last vertex {xN+i,yN+i) is assumed to 
be the same as the first one, so (xat+i, j/at+i) = { x \, yi ) 
and the polygon is closed. 

The body area of an individual fish was computed for 
every identified individual in every frame of a movie. All 
edge pixels forming a shape of an object (fish) were iden¬ 
tified based on the preset weighted intensity threshold. 
The resulting shape was a polygon with N vertices corre¬ 
sponding to a number of edge pixels with intensity values 
above the threshold. The area of this complex shape was 
computed with the same method as for the surface area 
of a group. 





FIG. S6. Experimental speed distributions for: (a) three sizes 
of fish in groups of 40 individuals, (b) small (1=7.5 mm), (c) 
medium (1=13 mm) and (d) large (1=23 mm) fish at variable 
group size. 


Additional motion statistics 


Figure IShl a) shows speed distributions for fish of three 
average body lengths at fixed group size (N=40). Small 
fish (1=7.5 mm) most of the time have speeds below 
100 mm/s and very often within a range of 0-20 mm/s. 
Very rarely small fish have speed above 200 mm/s. For 
medium-size individuals the distribution is much wider 
and has a peak at V = 125 mm/s. It spans up to V ss 250 
mm/s and has another minimum at 0-20 mm/s. For large 
fish the number of events when the individuals are sta¬ 
tionary or barely move decreases further and the peak 
is observed at V ~ 200. Figures ISOl bl-fd) show speed 
distributions for small, medium and large fish and vari¬ 
ous group sizes. All histograms for the same body size 
practically overlap. Therefore, while speed regime of fish 
depends strongly on body size it is not effected by the 
number of individuals in a group. 

Figure [ST] shows the average group area occupied by 
simulated fish as a function of number of individuals in 
a group. For all three sizes of simulated fish, the group 
area increases with group size, in agreement with experi¬ 
mental results (Fig. 3(a), main text). For all group sizes, 
big fish form densest groups. Conversely, groups of small 
fish occupy much larger area as compared to the other 
two size classes (medium-size and big fish). 
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FIG. S7. Average area of a simulated group of fish for three 
size classes: small (k = 1), medium-size [k = 1.73), and big 
(3.07) fish. Simulation parameters as stated in the main text. 
Methods section. 













